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As  one  of  the  key  components  of  proton-exchange  membrane  fuel  cells,  the  gas-diffusion  layer  (GDL) 
that  is  made  of  carbon  fibres  usually  exhibits  strong  structural  anisotropy.  Nevertheless,  the  GDL  has 
traditionally  been  simplified  as  a  homogeneous  porous  structure  in  modelling  the  transport  of  species 
through  the  GDL.  In  this  work,  a  coupled  electron  and  two-phase  mass  transport  model  for  anisotropic 
GDLs  is  developed.  The  effects  of  anisotropic  GDL  transport  properties  due  to  the  inherent  anisotropic 
carbon  fibres  and  caused  by  GDL  deformations  are  studied.  Results  indicate  that  the  inherent  structural 
anisotropy  of  the  GDL  significantly  influences  the  local  distribution  of  both  cathode  potential  and  current 
density.  Simulation  results  further  indicate  that  a  GDL  with  deformation  results  in  an  increase  in  the 
concentration  polarization  due  to  the  increased  mass-transfer  resistance  in  the  deformed  GDL.  On  the 
other  hand,  the  ohmic  polarization  is  found  to  be  smaller  in  the  deformed  GDL  as  the  result  of  the  decreased 
interfacial  contact  resistance  and  electronic  resistance  in  the  GDL.  This  result  implies  that  an  optimum 
deformation  needs  to  be  achieved  so  that  both  concentration  and  ohmic  losses  can  be  minimized. 

©  2008  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

The  proton-exchange  membrane  fuel  cell  (PEMFC)  is  regarded 
as  an  alternative  power  source  for  electric  vehicles  by  the  virtue  of 
its  high-energy  efficiency,  pollution-free  characteristics,  as  well  as 
the  simplicity  in  its  design  and  operation.  It  has  received  more  and 
more  attention  over  the  past  decade.  In  order  to  achieve  a  better 
performance,  it  is  essential  to  have  a  deep  understanding  of  the 
physicochemical  phenomena  occurring  in  each  component  of  the 
cell. 

As  one  of  the  key  components,  the  gas-diffusion  layer  (GDL), 
usually  made  of  highly  porous  material  such  as  carbon  fibre  papers 
or  cloths,  serves  to  facilitate  the  distribution  of  reactants  and 
removal  of  by-products,  and  to  provide  the  path  for  electrons  as 
well  as  heat  transport.  Since  GDLs  are  typically  made  of  fibrous 
porous  materials  they  exhibit  strong  structural  anisotropy  due  to 
the  special  constitutive  orientations  of  carbon  fibres.  This  intrin¬ 
sic  anisotropy  of  the  GDL  consequently  leads  to  different  transport 
properties  in  the  through-plane  and  in-plane  directions.  Among  the 
existing  PEMFC  models  reported  in  the  literature  however,  only  a 
few  [1-4]  consider  the  effect  of  the  inherent  anisotropy  of  the  GDL. 
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Zhou  and  Liu  [1]  numerically  studied  the  effect  of  the  electrical 
resistance  in  the  GDL  on  the  current  distribution  as  well  as  on  the 
cell  performance.  In  that  work,  the  anisotropic  characteristics  of 
the  GDL  were  evaluated  by  using  different  electrical  conductivities 
in  the  through-plane  and  in-plane  directions.  Pasaogullari  et  al.  [2] 
presented  a  two-dimensional  model  to  investigate  the  effects  of  the 
anisotropic  properties  on  the  coupled  heat  and  water  transport  in 
the  GDL.  Their  results  indicated  that  the  inherent  anisotropy  of  the 
GDL  plays  an  important  role  in  determining  the  temperature  distri¬ 
bution  in  the  GDL.  Recently,  Pharoah  et  al.  [3]  briefly  reviewed  the 
approaches  currently  employed  to  define  the  effective  transport 
coefficients  in  the  GDL.  They  also  developed  a  two-dimensional 
single-phase  model  to  examine  the  effect  of  the  inherent  anisotropy 
of  the  GDL  on  the  profile  of  the  local  current  density  and  the  cell  per¬ 
formance.  Although  both  the  isotropic  model  and  the  anisotropic 
models  gave  nearly  identical  polarization  curves  for  certain  sets  of 
parameters,  the  profiles  of  local  current  density  are  quite  different. 
More  recently,  a  simplified  method  was  also  developed  by  Meng 
[4]  for  solving  the  anisotropic  transport  problem  in  a  PEMFC. 

In  addition  to  the  inherent  anisotropy,  the  GDL  usually  under¬ 
goes  deformation  when  all  the  cell  components  are  assembled.  This 
deformation  causes  significant  changes  in  the  physical  properties 
of  the  GDL  and  their  distributions.  Since  the  GDL  tends  to  have  a 
larger  deformation  under  the  rib  region  than  under  the  channel 
region,  inevitably,  the  local  properties  of  the  GDL  show  strongly 
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Nomenclature 

Ay  specific  area  (m-1 ) 

C  molar  concentration  (mol  m-3 ) 

D  diffusivity  (m2  s-1 ) 

F  Faraday  constant  (96  478  C  mol-1) 

i’o  exchange  current  density  (A  m-3 ) 

i+  proton  current  vector  (A  m-2 ) 

/  current  density  (A  m-2 ) 

jc  cathodic  current  density  (A  m-3 ) 

j p  source  term  in  proton  conservation  equation 

(Am-3) 

js  source  term  in  electron  conservation  equation 

(Am-3) 

kc  condensation  rate  (s-1 ) 

ke  evaporation  rate  (atm-1  s-1 ) 

/<H  Henry’s  law  constant 

kr  relative  permeability 

I<  permeability  of  porous  material  (m2 ) 

L  thickness  of  porous  layers  (m) 

m  source  term  in  mass  conservation  equation 

(kgm-3  s-1) 

mpt  Pt  loading  in  cathode  catalyst  layer  (kg  m-2 ) 

M  molecular  weight  (kg  mol-1 ) 

nd  electro-osmotic  drag  coefficient 

N  mol  flux  (mol  m-2  s-1) 

pc  capillary  pressure  (Pa) 

pg  gas-phase  pressure  (Pa) 

Pi  liquid-phase  pressure  (Pa) 

R  gas  constant  or  radius  (J  mol-1  K-1 )  or  (m) 

R  source  term  in  species  conservation  equation 

(molm-3  s-1) 

R  interfacial  species  transfer  rate  (mol  m-3  s-1 ) 

Rc unt  interfacial  contact  resistance  ( £2  m2 ) 

s  liquid  saturation 

T  temperature  (I<) 

Uq  thermodynamic  equilibrium  voltage  of  ORR  (V) 

Vc  cathode  voltage  (V) 

x  coordinate  (m),  or  mole  fraction  in  liquid  solution 

(mol  mol-1) 

y  coordinate  (m),  or  mole  fraction  in  gas  mixture 

(mol  mol-1) 

Greek  letters 

ac  charge-transfer  coefficient  at  cathode 

8  thickness  of  Nation  film  (m) 

s  porosity  of  porous  medium 

r\  overpotential  (V) 

0C  contact  angle  (°) 

/x  viscosity  (kg  m-1  s-1 ) 

p  density  (kgm-3) 

a  interfacial  tension  (N  m-1 ) 

cre  ionic  conductivity  of  membrane  (£2-1  m-1 ) 

gs  electrical  conductivity  (£2-1  m-1 ) 

Subscripts 
c  cathode 

cl  catalyst  layer 

g  gas  phase 

gdl  gas  diffusion  layer 

in  inlet  condition 

in-p  in-plane  direction 

1  liquid  phase 


mem  membrane 

N  Nation  phase 

th-p  through-plane  direction 
W  water 

WV  water  vapour 

Superscripts 

eff  effective  value 

ref  reference  value 

sat  saturated  value 


inhomogeneous  distributions  along  the  in-plane  direction.  In  addi¬ 
tion,  the  interfacial  contact  resistance  is  also  changed,  which  is 
known  as  one  of  the  important  parameters  that  affects  the  voltage 
loss  in  a  PEMFC.  Experimental  studies  [5-8]  have  been  conducted 
to  examine  the  effect  of  the  deformation  on  the  physical  properties 
of  the  GDL,  such  as  porosity,  electrical/thermal  conductivity  and 
permeability,  as  well  as  interfacial  contact  resistance.  The  effect 
of  the  GDL  deformation  on  the  performance  of  a  PEMFC  has  also 
been  experimentally  investigated  [9,10].  To  date,  few  models  have 
been  reported  to  examine  the  effect  of  GDL  deformation.  Zhou  et  al. 
[11,12]  numerically  studied  the  effect  of  the  GDL  deformation  on  the 
performance  of  a  PEMFC.  Although  the  influence  of  GDL  deforma¬ 
tion  on  the  contact  resistance  and  porosity  distribution  were  taken 
into  account  in  their  models,  the  impact  on  electrical  resistance  in 
the  GDL  was  not  considered.  Furthermore,  the  anisotropic  nature  of 
the  GDL  was  neglected  in  their  model.  Hottinen  et  al.  [13]  numer¬ 
ically  studied  the  influence  of  GDL  deformation  on  the  profiles  of 
local  current  density  by  using  a  single-phase  model.  Although  the 
anisotropic  properties  of  the  GDL  such  as  electronic  conductivity 
and  permeability  were  taken  into  account,  the  anisotropic  diffu¬ 
sion  coefficients  of  oxygen  transport  in  the  GDL  was  not  indicated. 
More  recently,  Hottinen  et  al.  [14]  studied  the  effect  of  the  GDL 
deformation  on  the  temperature  distribution  in  a  PEMFC. 

The  above  review  of  the  literature  indicates  that  there  is  still 
a  lack  of  a  sound  mathematical  model  that  can  take  into  account 
the  influences  of  both  the  intrinsic  GDL  anisotropy  associated  with 
carbon  fibres  and  the  structural  GDL  anisotropy  caused  by  inhomo¬ 
geneous  GDL  deformation.  The  objective  of  the  work  is  to  develop  a 
2D  model  that  can  simulate  coupled  electron  and  two-phase  mass 
transport  in  a  PEMFC  cathode.  With  this  model,  we  numerically 
have  studied  the  effects  of  anisotropic  GDL  properties  associated 
with  carbon  fibres  and  caused  by  GDL  deformation  on  the  distribu¬ 
tions  of  reactant,  potential,  overpotential  and  the  current  density 
as  well  as  on  the  cell  performance. 


2.  Formulation 


Consider  the  physical  domain  sketched  in  Fig.  1,  which  consists 
of  a  cathode  GDL,  a  cathode  catalyst  layer  (CL)  and  a  polymer  elec¬ 
trolyte  membrane  (PEM).  Note  that  since  both  channel  width  and 
rib  width  in  a  parallel  flow-field  are  symmetrical  with  respect  to 
their  middle  points,  only  a  half-rib  width  and  a  half-channel  width 
need  to  be  considered  to  save  the  computing  time.  The  mass  trans¬ 
port  of  gas  oxygen  and  water  vapour  in  the  GDL  are  governed  by 
the  general  convection-diffusion  transport  equation,  i.e. 


d  .  d  d 

^(ugQ,g)+^(ygQ,g)-^ 


dx 


dQ,g 

dx 


d_ 

dy 


dQ,gA 

) 


=  Rh  i:  02,WV 


(1) 
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Fig.  1.  Modelling  domain:  (a)  geometry  of  the  GDL  without  deformation;  (b)  geometry  of  the  GDL  with  deformation.  Roman  numerals  refer  to  boundaries  of  modelled 
geometries. 


where  Df^  is  the  effective  diffusion  coefficient  of  species  in  a  porous 
medium.  For  a  porous  structure  that  is  made  of  random  fibrous 
materials,  like  GDLs,  Df^  is  given  by  [15,16] 

Dfg  =  D°/(g)g(s)  -  D°g£(i  1 0.11  )g(1  ~s)fi  (2) 


phases  is  embodied  in  the  source  terms  on  the  right-hand  sides  of 
Eqs.  (4)  and  (5). 

With  respect  to  the  dissolved  water  in  the  membrane,  water 
transfer  is  driven  by  the  concentration  gradient  and  electro-osmotic 
drag.  The  general  conservation  of  dissolved  water  in  the  membrane 
can  be  expressed  by 


where  the  constant  a  is  0.521  and  0.785  for  the  in-plane  and 
through-plane  diffusions,  respectively.  The  constant  is  set  to  be 
1.5  for  both  types  of  diffusion.  In  the  CL,  the  effective  diffusion 
coefficient  of  a  given  species  can  be  obtained  from  the  Bruggeman 
correlation: 


,0  c1.5 


DZ=DU£ 


(1  -sf 


(3) 


In  Eq.  (1),  ug  and  vg  represent  the  superficial  velocities  of  gas 
mixture  in  a  porous  medium,  that  can  be  calculated  by  Darcy’s  law 
if  the  pressure  field  is  available.  The  gas  pressure  can  be  obtained 
by  solving 


d_ 

dx 


Kkrg  dpg  \  9  (_  Kkrg  dpg 

g  Ug  dx  J  +  dy\  g  Mg  dy 


=  mg 


(4) 


The  liquid  saturation  distribution  in  the  cathode  porous  region, 
s,  can  be  modelled  by 


d_  [/  Pi  Mri  \  (  dpc  3s,  9pg  \ 
dx  v  /Xi  )  y  ds  dx  dx  J 


(8) 


where  the  electro-osmotic  drag  (EOD)  coefficient,  nd,  and  the  dif- 
fusivity  of  water  in  the  membrane,  Dw, n»  are  both  functions  of  the 
water  content  in  the  membrane. 

In  addition  to  the  two-phase  mass  transport  processes  in  the 
porous  cathode,  the  transport  of  electrons  through  the  catalyst  layer 
and  gas-diffusion  layers  and  the  transport  of  protons  through  the 
membrane  phase  must  also  be  considered.  The  governing  equa¬ 
tions  that  describe  the  transport  of  electrons  and  protons  in  the 
membrane  phase  are,  respectively,  given  by 


(9) 


and 


d 


dpc  9si  9pg 
ds  dy  +  dy 


=  m. 


(5) 


where  pc  is  the  capillary  pressure,  which  represents  the  difference 
between  the  gas  pressure  and  liquid  pressure  as  follows: 


/  £  \  ®  ^ 

Pc=Pg-Pi  =0  COS0c(  J(s) 


(6) 


In  order  to  account  for  the  effect  of  water  evaporation  and  con¬ 
densation,  the  interfacial  transfer  rate  of  water  between  the  liquid 
and  gas  phase  is  given  by  [17] 


fee  (ywvPg  -  Pwv  ),  ywvPg  <  Pw 

M[H20 

fee ~Rs^ywv  (ywvPg  -  Pwv),  ywvPg  >  Pw 


(7) 


(10) 


where  crs  is  the  electrical  conductivity  of  the  electron-conducting 
phase,  while  cre  is  the  proton  conductivity  in  the  membrane,  which 
is  known  as  to  be  a  function  of  water  content  in  the  membrane. 

Up  to  this  point,  we  have  presented  all  the  governing  equations 
that  describe  two-phase  mass  transport  and  charge  transport  in  the 
PEMFC  cathode.  To  make  the  above  governing  equations  closed, 
some  constitutive  correlations  and  definitions  are  needed.  These 
include  capillary  pressure,  relative  permeability  for  both  phases, 
sources  terms,  and  so  on.  All  these  correlations  and  associated 
nomenclatures  are  listed  in  Table  1. 


2.1.  Boundary  conditions 


where  andywv  denote  the  saturation  pressure  of  water  vapour 
and  the  mole  fraction  of  water  vapour  in  the  cathode  gas  phase, 
respectively.  Note  that  the  interfacial  transfer  of  water  between  the 


As  indicated  in  Fig.  1 ,  there  are  seven  boundaries  specified  in  the 
modelling  domain  as  marked  with  Roman  numerals.  The  conditions 
at  each  boundary  are  described  below. 
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Table  1 

Constitutive  correlations  and  definitions 


Parameters 


Expressions 


Capillary  pressure 

Relative  permeabilities 

General  generation  rate  of  mass  in  liquid 
phase 

General  generation  rate  of  mass  in  gas 
phase 

Mole  generation  rate  of  species 

General  generation  rate  of  charges 

Thickness  of  the  GDL  under 
inhomogeneous  compression  [5] 

Porosity  of  the  GDL  under  inhomogeneous 
compression  [5] 

Electrical  conductivity  of  the  GDL  [5] 

Interfacial  contact  resistance  [5] 


Pc=Pg-Pi=cr cos 0c(£/K)o'5;(s),i(s)  =  {  _  2.120s2  + 

/<rl  =  s3  liquid;  krg  =  (1  -  s)3 


90°  <  6c  <  180° 


gas 


Mh20^w 
Jh2o  jp  ~ 

Mh20^w 


GDL 


h?i  =  <  jc  ..  Nw,  cross  ~ 

'  MH2o  2p  +  mh20  — 2^ - mh2o4  CL 


-Mq2  ^  +Mh2o4  CL 


{ 

|_|L’  ^wv,c  =  {| 

_  /  0  MEM  /  0 

~\-jc  cl  ’  Js  “  \  Jc 


Rq2  = 


J  p 


Rw  GDL 
Rw  CL 


GDL 

CL 


{Lcoi 

19. 


30314  log((y  -  0.0005)  x  106  +  1 )  x  1CT6  +  Lcomp 


,  Lcomp :  thickness  of  GDL  under  the  rib  (LCOmp  =  250  p.m) 


^(y )  =  ^0  Lgdi-L^in  ’  Lmin  =  ^ _  S° 

<xin_p(y)  =  6896  -  1.159  X  107  L(y),  <xth_p(y)  =  3285  -  8.385  x  106  L(y) 
l?cont(y)  =  5.83  x  10-10  exp[2.06  x  104  L(y)] 


Boundary  I.  This  boundary  represents  the  inlet  of  the  reactant 
supply  in  the  cathode,  at  which  the  concentration  of  gas  species, 
gas  phase  pressure,  and  liquid  water  saturation  are  all  specified 
to  be  inlet  conditions.  With  respect  to  the  electrical  potential,  the 
gradient  is  set  to  be  zero  at  the  direction  normal  to  this  boundary: 

Q)2,g  =  Q)2,in>  CWV,g  —  Q/w,in?  Pg  —  Pg,in>  s  —  ^im 

I  normal  to  I  —  0  (H) 


where  aw  represents  the  water  vapour  activity  (aw  = 

*vaporPg,a/  Pvapor)- 

Boundary  VI.  This  boundary  represents  the  interface  between 
the  GDL  and  the  CL.  The  electrical  current  and  the  fluxes  of  oxy¬ 
gen,  vapour  and  liquid  water  are  continuous  at  this  interface  to 
satisfy  the  general  mass/species/charge  balance.  Since  there  is  no 
ionic  phase  in  the  GDL,  the  flux  of  protons  is  set  to  be  zero  at  this 


Boundary  II.  This  boundary  is  the  interface  between  the  GDL  and 
the  rib  collector,  which  is  an  impermeable  wall  to  mass/species 
transfer.  Accordingly,  the  flux  of  the  gaseous  species  and  liquid 
species  in  the  x  direction  (normal  to  boundary  II)  are  specified  to  be 
zero.  At  this  interface,  the  electrical  potential  is  set  to  be  the  applied 
potential,  i.e. 


dQ>2.g 

dx 


=  0, 


9cwv,g  0  dPg  _n 
dx  '  dx 


|  =  °,  4>s  =  vc 


(12) 


Boundaries  III  and  IV.  These  two  boundaries  are  symmetrical  to 
the  middle  point  of  rib  width  and  the  middle  point  of  channel  width, 
respectively.  Hence,  the  gradients  of  all  the  variables  in  x  direction 
are  set  to  zero: 


Electrolyte  and 
Pt/Carbon  particles 

Electrolyte  film  covering 
the  agglomerate 

Void  space  occupied  by 
gas  and  water 


9Co2.g  =  Q  9Cwv,g  =  Q 
dy  ’  dy 


9pg=0  9cw,n  Q 

dy  ’  dy  ’  dy 

(13) 


Boundary  V.  This  boundary  represents  the  interface  between  the 
membrane  and  the  anode.  At  this  boundary,  it  is  assumed  that  the 
water  concentration  is  only  a  function  of  the  water  vapour  activity 
in  the  anode  gas  stream  and  an  equilibrium  condition  is  reached 
between  this  water  vapour  and  the  dissolved  water  in  the  mem¬ 
brane.  As  for  the  ionic  potential,  the  latter  is  set  to  zero  as  a  reference 
at  this  boundary: 

Cw,N  =  CNa,  </>e  =  0  (14) 

The  equilibrium  water  concentration  when  the  membrane  is 
only  in  contact  with  water  vapour  [18]  is  obtained  from 

C^N  =  Cf( 0.043  +  17.81aw  -  39.85a^,  +  36.0o£,)  (15) 


r=0  r=Ragg 

Fig.  2.  Schematic  of  agglomerate  structure  in  cathode  catalyst  layer. 


Table  2 

Correction  factors  and  correlations  in  thin  film-agglomerate  model 


Parameters,  symbols 

Expressions 

Correction  factor,  ^ 

Correction  factor,  £2 

^  1+[((5n/(<5n  +^agg ) )  1  (Pq^  agg /^°2 > N ^ ^a§g  coth( 0^agg )~1] 

£2  _  (0J?  )2  I^agg  COth(0Ragg)  1] 

Correction  factor,  £3  [21] 

=  1  -  S 

The  Thiele  modulus,  0Ra gg 

Effective  specific  area  in 
agglomerate,  Aagg 

Effective  diffusion 
coefficient  of  oxygen, 

03,agg 

0Kagg  -  \  nFDeff  Cref  eXP  (  2RT  ^c)  RaSg 

y  02,agg  02 

A  —  A 

agg  £s+£N,agg(£s/(l-£N,agg)) 

neff  _  XIpt-5 

L)2,agg  °2>N£N,agg 

Table  3 

Cell  geometric  dimensions  and  operating  conditions 
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Parameters 

Symbols 

Value 

Unit 

Thicknesses  of  functional  layers 

Tgdl  -  tc\ ,  Lmem 

0.38,  0.02,  0.05 

mm 

Width  of  channel,  rib 

Wc,  Wr 

1.0,  1.0 

mm 

Operation  temperature 

T 

348.15 

K 

Inlet  pressure  of  cathode 

P  g.in 

1.013  x  105 

Pa 

Inlet  oxygen  concentration 

Co2,in 

7.35 

mol  m3 

Inlet  water  vapour  concentration 

Cwv.in 

0 

mol  m3 

Inlet  liquid  saturation  of  cathode 

Sc,in 

0 

- 

Relative  humidity  of  anode 

RH 

100% 

- 

boundary: 


22.  Electrochemical  kinetics 


Nt |-  =  Nj|+,  i :  02,WV,W, 


(16) 


Boundary  VII.  This  boundary  represents  the  interface  between 
the  CL  and  the  PEM.  At  this  interface,  it  is  simply  assumed  that  the 
water  concentration  in  the  membrane  is  in  equilibrium  with  the 
water  in  the  catalyst  layer  and  varies  linearly  as  a  function  of  liquid 
saturation  in  the  catalyst  layer  [18]: 


C$  N  =  14Cf  +  2.8Cfs  if  s  >  0  (17) 

Due  to  the  nature  of  the  membrane,  the  fluxes  of  oxygen,  water 
vapour  and  electronic  current  in  the  x  direction  are  set  to  zero  at 
this  interface,  where  the  flux  of  liquid  water  and  the  ionic  current 
are  continuous,  i.e. 


dC02,S  0  dcwv,g 

dx  ’  dx 


9pg 

dx 


=  0, 


Nw\-  =  Nw|+, 


r  rec 1 

Cw,N  =  Cw  N,  ae 


■■  <7e 


90e 

dx 


(18) 


With  respect  to  the  electrochemical  oxygen  reduction  reaction 
(ORR)  at  the  cathode  catalyst  layer,  the  transfer  current  based  on 
first-order  Tafel  kinetics  can  be  written  as 

ic=Avi0(^f  exp(!^c)  (19) 

°2 

where  Co2  is  the  surface  concentration  of  oxygen  at  the  active  reac¬ 
tion  sites  (i.e.  three-phase  interfaces),  and  rjc  represents  the  cathode 
overpotential  which  is  calculated  by 

Be  —  Uo  ~  0s  +  0e  (20) 

It  has  been  well  documented  that  the  active  reaction  sites 
are  usually  surrounded  either  by  electrolyte  films  or  by  liquid 
water  films  when  a  two-phase  flow  exists.  Oxygen  in  gas  pores 
needs  to  permeate  through  these  media  to  the  active  reaction 
sites.  It  is  thus  necessary  to  consider  the  additional  mass-transfer 
resistance  of  oxygen  from  gas  pores  to  the  actual  reaction  sites 
because  the  diffusivity  of  oxygen  in  Nation  or  liquid  water  is 
quite  low.  Conventionally,  a  spherical  thin-film  agglomerate  model 


Table  4 

Physicochemical  properties 


Parameters,  symbols 

Value 

Unit 

Reference 

Porosity  of  GDL,  £gdj 

0.7 

- 

[19] 

Porosity  of  membrane,  £mem 

0.3 

- 

[19] 

Nafion  volume  fraction  in  the  agglomerate,  £agg 

0.5 

- 

Volume  fraction  of  solid  in  CCL,  £s 

0.3 

- 

Radius  of  the  agglomerate,  Ragg 

0.6 

pan 

Thickness  of  Nafion  film  covering  the  agglomerate,  5n 

20 

nm 

Permeability  of  GDL,  /Cth_p,  I<in_p 

8.69  xlO-12,  3.0  xlO-12 

m2,  m2 

[2] 

Permeability  of  CL,  I< 

1.0  x  10“14 

m2 

[19] 

Permeability  of  MEM,  K 

2.0  x  10“18 

m2 

[19] 

Electric  conductivity  of  the  GDL,  crs  thrp,  crs,inrp 

100,  2500 

C2-1  m_1 

[5] 

Electric  conductivity  of  CL,  os 

90 

C2_1  m_1 

[1] 

Diffusivity  of  02  in  gas,  D02,g 

1.775  x  10-5(T/273.15)1'823 

m2  s_1 

[21] 

Diffusivity  of  02  in  water,  Dq2\ 

3.032  x  10“9 

m2  s-1 

[22] 

Diffusivity  of  02  in  Nafion,  Dq2,n 

1.844  x  10-10 

m2  s-1 

[22] 

Diffusivity  of  water  vapour,  Dwv,g 

2.56  x  10-5(T/307.15)2'334 

m2  s-1 

[21] 

Diffusivity  of  water  in  the  membrane,  DWN 

4.17  x  10~8A,(161  e“*  + 1 )  e~2436lT 

m2  s_1 

[24] 

Proton  conductivity,  cre 

(0.5139A.  -  0.326)  x  exp[  1268.0(1.0/303  -  1.0/7)] 

£2-1  m_1 

[18] 

Water  content  in  membrane,  1 

Cw,n/Q 

- 

[18] 

Fixed  charge  sites  concentration  in  membrane,  Cf 

1200 

mol  m-3 

[18] 

Coefficient  of  electro-osmotic  drag  of  water,  nd 

(2.5/22)1 

- 

[18] 

Viscosity  of  gas,  liquid  phase,  /zg,  /q 

2.03  xlO-5,  4.05  xlO-4 

kgm-1  s_1 

[15] 

Evaporation  rate  constant  for  water,  ke 

1 

atm-1  s_1 

[23] 

Condensation  rate  constant  for  water,  kc 

100 

s-1 

[23] 

Henry  law  constant  for  oxygen,  kH,o2 

e(-666/T+14.1 )  !pjT 

- 

[22] 

The  saturation  pressure  of  water  vapour,  log10  p ^ 

-2.1794  +  0.02953(7-  273)  -  9.1837  x  10-5(7-  273 )2  + 1.4454  x  10-y(7-  273 )3 

atm 

[23] 

Thermodynamic  voltage  of  ORR,  Uo2 

1.21 

V 

Pt  loading  in  CCL,  mPt 

0.4 

mg  cm-2 

Specific  area  in  CCL,  Av 

A0mPtlLCi 

ITT1 

[22] 

Exchange  current  density  of  ORR,  i0 

10(3.507—4001  IT) 

A  cm-2 

[22] 

Transfer  coefficient,  ac 

1.0 

- 

[22] 

Reference  concentration  of  oxygen,  CJ®f 

1.2 

mol  m-3 

[22] 

770 


W.W.  Yang  et  al.  /  Journal  of  Power  Sources  185  (2008)  765-775 


[19,20]  is  employed  to  capture  the  influence  of  the  microstructure 
of  catalyst  layer,  as  sketched  in  Fig.  2.  With  further  taking  into 
account  the  additional  mass-transfer  resistance  and  the  interfa¬ 
cial  mass  transfer  resistance,  the  kinetics  of  the  ORR  is  expressed 
as 

ic=Av/o  exp(^>;c)?i?2?3  (21) 

where  /<h»  §i  ,  £2  and  §3  denote  the  Henry  factor  to  capture  the  effect 
of  the  dissolving  process,  the  correction  factor  in  view  of  the  trans¬ 
fer  resistance  in  the  external  Nafion-film,  the  correction  factor  in 
view  of  the  transport  resistance  in  the  agglomerate,  and  the  cor¬ 
rection  factor  to  capture  the  effect  of  the  liquid  water  surrounding 
agglomerate,  respectively.  The  detailed  expressions  of  the  above 
factors  and  correlations  are  listed  in  Table  2. 


.300 

.296 
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CD 
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CD 
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.288 
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> 

O 
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.284 

lb 

3.  Results  and  discussion 

The  governing  equations  for  the  cell  geometric  dimensions  and 
operating  parameters  listed  in  Table  3  subjected  to  the  boundary 
conditions,  along  with  the  constitutive  relations  shown  in  Table  1 
and  the  physicochemical  properties  listed  in  Table  4,  are  solved 
numerically  using  a  self-written  code,  which  was  developed  based 
on  the  SIMPLE  algorithm  with  the  finite-volume-method  [19,20]. 
Once  the  profile  of  the  ionic  potential  was  obtained,  the  local  cur¬ 
rent  density  was  calculated  based  on  the  flux  of  protons  across 
the  CL/PEM  interface,  and  the  average  current  density  was  then 
estimated  by 

_  /oWc/2+Wr/2^e(30e/9x)|CL/MEMdy  _  JJCLjcdxdy 
cell_  Wr/2  +  Wc/2  ~  Wx/2  +  Wc/2  1  j 

In  this  work,  simulations  were  performed  for  three  different 
cases.  Case  I:  the  GDL  is  treated  as  an  isotropic  porous  medium 
and  transport  coefficients  in  the  GDL  are  based  on  its  through-plane 
properties,  as  sketched  in  Fig.  1  (a).  Case  II:  the  GDL  as  an  anisotropic 
structure  without  deformation,  also  as  sketched  in  Fig.  1  (a).  Case  III : 
the  GDL  as  an  anisotropic  structure  with  deformation. 

3.1.  Isotropic  GDL 

The  profiles  of  oxygen  concentration,  local  current  density  and 
potentials  in  the  CL  along  they  direction  at  a  given  cathode  potential 
of  0.4  V  are  shown  in  Fig.  3.  Due  to  the  rib  coverage  effect,  oxygen 
concentration  is  lower  in  the  region  under  the  rib  than  that  in  the 
region  under  the  channel.  Higher  oxygen  concentration  will  natu¬ 
rally  result  in  higher  cell  performance.  As  shown  in  Fig.  3(a),  as  a 
whole,  the  current  density  is  greater  in  the  region  under  the  channel 
than  that  in  the  region  under  the  rib.  It  is  interesting  to  note,  how¬ 
ever  that  the  highest  current  density  occurs  in  the  location  near  the 
interface  between  the  rib  and  channel  (i.e.,  x  =  0.55  mm),  although 
the  highest  oxygen  concentration  occurs  in  the  location  under  the 
middle  of  the  channel  (i.e.,  x  =  1.0  mm).  This  can  be  explained  by  the 
fact  that  the  local  current  density  is  influenced  by  both  the  local 
oxygen  concentration  and  the  cathode  overpotential  in  the  CL.  As 
shown  in  Fig.  3(b),  the  cathode  overpotential  decreases  gradually 
from  the  region  under  the  rib  to  the  region  under  the  channel,  as 
opposed  to  the  variation  of  the  oxygen  concentration  in  the  CL. 
Both  variations  of  the  oxygen  concentration  and  cathode  overpo¬ 
tential  affect  the  current  density  distribution.  The  ultimate  profile 
of  current  density  in  the  CL  is  determined  by  the  combined  effects  of 
oxygen  concentration  and  overpotential  changes.  The  profiles  of  the 
electrical  potential  and  ionic  potential  in  the  CL  are  also  displayed  in 
Fig.  3(b).  The  electric  potential  gradient  is  to  drive  the  electric  cur¬ 
rent  from  all  the  regions  inside  the  GDL  and  CL  to  the  rib  collector. 
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Fig.  3.  Variations  in  oxygen  concentration  and  local  current  density  (a),  electrical 
potential,  ionic  potential  and  overpotential  in  the  CL  along  y  direction  (b)  at  a  given 
cathode  potential  of  0.4  V. 


Hence,  in  the  regions  away  from  the  rib  collector,  larger  potential 
gradients  are  required.  This  is  why  the  electrical  potential  is  always 
higher  in  the  region  under  the  middle  of  the  channel,  as  shown  in 
Fig.  3(b).  Notice  that  the  difference  between  the  electric  potential 
in  the  CL  and  the  cathode  voltage  (i.e.,  the  electrical  potential  at  the 
ribvertical  live  GDL  interface)  represents  the  ohmic  loss  due  to  the 
resistance  of  the  GDL.  With  respect  to  the  ionic  potential  in  the  CL,  it 
is  more  uniform  along  the  y  direction  compared  with  the  distribu¬ 
tion  of  the  electric  potential.  With  setting  a  reference  value  of  zero 
at  the  vertical  live  PEM/CL  interface,  the  ionic  potential  is  always 
negative  in  the  cathode  catalyst  layer.  The  absolute  value  of  ionic 
potential  herein  denotes  the  ohmic  loss  across  the  membrane  as  a 
result  of  the  proton  transfer  resistance. 

The  variation  in  potential  along  the  x  direction  is  also  displayed 
in  Fig.  4.  A  relatively  large  voltage  loss,  approximate  150  mV,  across 
the  membrane  can  be  observed  compared  with  the  voltage  loss, 
about  50  mV,  across  the  diffusion  layer.  This  is  because  the  ohmic 
resistance  in  the  membrane  is  much  higher  than  that  in  the  GDL. 

3.2.  Anisotropic  GDL  without  deformation 

The  GDL  exhibits  significant  anisotropy,  leading  to  differ¬ 
ent  transport  properties  along  the  in-plane  and  through-plane 
directions.  To  evaluate  the  influence  of  the  anisotropic  electrical 
conductivity  of  the  GDL  on  the  distributions  of  electrical  potential 
and  local  current  density,  simulations  were  conducted  by  keeping 
the  same  through-plane  electrical  conductivity  of  the  GDL  but  vary¬ 
ing  the  in-plane  electrical  conductivity.  The  profiles  of  the  electrical 
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Through-plane  direction  /  mm 

Fig.  4.  Variations  in  electric  potential,  ionic  potential,  and  overpotential  along  x 
direction  for  given  cathode  potential  of  0.4  V. 


Fig.  6.  Variations  in  local  current  density  along  y  direction  for  various  in-plane 
electric  conductivities  of  GDL  at  given  cathode  potential  of  0.4  V. 


potential  in  the  CL  along  the  y  direction  for  various  in-plane  electri¬ 
cal  conductivities  are  presented  in  Fig.  5.  As  the  in-plane  electrical 
conductivity  is  increased,  the  electron  transport  along  the  in-plane 
direction  is  enhanced,  thereby  facilitating  the  transport  of  the  elec¬ 
trical  current  from  the  region  under  the  channel  to  the  rib  collector. 
As  shown  in  Fig.  5,  clearly,  the  electrical  potential  in  the  CL  decreases 
significantly  as  the  in-plane  electric  conductivity  is  increased,  espe¬ 
cially  that  in  the  region  under  the  channel.  In  turn,  the  decrease  in 
the  electrical  potential  in  the  CL  results  in  an  increase  in  the  rate  of 
current  generation  in  the  CL.  The  profiles  of  the  local  current  density 
along  the  y  direction  for  various  in-plane  electrical  conductivities 
are  shown  in  Fig.  6.  As  a  whole,  higher  in-plane  electric  conduc¬ 
tivity  leads  to  higher  current  density.  The  increase  in  the  current 
density  is  more  pronounced  in  the  region  under  the  channel  than 
that  in  the  region  under  the  rib,  as  the  electrical  potential  is  more 
significantly  decreased  in  the  region  under  the  channel  than  that 
in  the  region  under  the  rib. 

The  effective  diffusion  coefficient  of  species  in  a  fibrous  porous 
structure  also  exhibits  anisotropic  behaviour.  For  example,  in  a 
fibrous  porous  medium  with  a  porosity  of  0.7,  the  effective  fac¬ 
tor/^)  in  Eq.  (2)  is  0.51  in  the  through-plane  direction,  but  0.57  in 
the  in-plane  direction.  To  examine  the  influence  of  the  anisotropic 
diffusion  coefficients  of  oxygen  in  the  GDL,  an  isotropic  GDL  hav¬ 


ing  the  through-plane  effective  diffusion  coefficient  is  compared 
with  an  anisotropic  GDL  having  different  diffusion  coefficients  in 
the  through-plane  and  in-plane  directions.  Fig.  7  shows  the  pro¬ 
files  of  local  oxygen  concentration  and  local  current  density  along 
the  y  direction  at  a  given  cathode  potential  of  0.4  V.  Clearly,  the 
oxygen  concentration  in  the  CL  predicted  by  the  anisotropic  trans¬ 
port  model  is  higher  than  that  predicted  by  the  isotropic  transport 
model,  especially  in  the  region  under  the  rib.  This  is  primarily  due  to 
the  fact  that  the  higher  in-plane  diffusion  coefficient  benefits  the 
transport  of  oxygen  in  the  y  direction.  Consequently,  the  current 
density  is  also  higher  when  the  anisotropic  transport  of  oxygen  in 
the  GDL  is  taken  into  account,  as  shown  in  Fig.  7. 

The  fuel  cell  performance  predicted,  by  the  isotropic  and 
anisotropic  transport  models  is  compared  in  Fig.  8.  At  low-current 
densities  (<0.3  A  cm-2),  the  cell  performance  predicted  by  the  two 
models  is  almost  the  same,  as  both  the  ohmic  and  the  concentration 
polarization  at  low-current  densities  are  negligible.  With  increase 
in  cell  current  density,  a  voltage  difference  between  the  two  curves 
can  be  observed  and  becomes  significant  at  large  current  densities. 
For  example,  at  a  current  density  of  1.0  A  cm-2,  the  voltage  pre¬ 
dicted  by  the  anisotropic  transport  model  is  about  25  mV  higher 
than  that  predicted  by  the  isotropic  transport  model.  The  data  in 
Fig.  8  also  show  that  the  maximum  power  density  predicted  by 
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Fig.  5.  Variations  in  electrical  potential  in  CL  along  y  direction  for  various  in-plane  Fig.  7.  Variations  in  oxygen  concentration  and  current  density  in  CL  alongy  direction 

electric  conductivities  of  GDL  at  given  cathode  potential  of  0.4  V.  for  both  isotropic  transport  and  anisotropic  GDLs. 
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Fig.  8.  Comparison  in  fuel  cell  performance  simulated  with  isotropic  transport 
model  (case  A)  and  anisotropic  transport  model  (case  B). 


the  anisotropic  model  is  about  9%  higher  than  that  by  the  isotropic 
model. 

3.3.  Anisotropic  GDL  with  deformation 

Deformation  of  the  GDL  is  usually  encountered  in  a  fuel  cell 
when  all  the  cell  components  are  assembled.  Due  to  the  deforma¬ 
tion,  the  thickness  and  the  physical  properties  of  the  GDL  and  the 
interfacial  contact  resistance  are  all  changed.  For  the  purpose  of 
modelling,  the  set  of  correlations  given  in  Ref.  [  5  ],  as  listed  in  Table  1 , 
is  simply  employed  to  capture  the  effect  of  the  GDL  deformation. 
It  is  particularly  worth  mentioning  that  the  GDL  deformation  is 
not  homogeneous.  The  parts  of  the  GDL  situated  under  the  rib  are 
more  significantly  compressed  than  those  under  the  channel.  Such 
inhomogeneous  deformation  may  cause  significant  variation  in  the 


local  physical  properties  of  the  GDL  along  the  y  direction,  thereby 
influencing  local  current  densities. 

To  examine  the  effects  of  GDL  deformation  on  the  reactant, 
potentials,  and  current  distributions,  as  well  as  on  the  overall  cell 
performance,  the  following  three  cases  are  compared:  (a)  isotropic 
GDL  without  deformation,  (b)  anisotropic  GDL  without  deforma¬ 
tion  and  (c)  anisotropic  GDL  with  deformation. 

'S'  The  distribution  of  oxygen  concentration  in  the  GDL  and  CL  for 

the  three  cases  are  given  in  Fig.  9.  The  simulations  are  performed 
at  the  same  current  density  of  1.1  A  cm-2.  A  comparison  between 
Fig.  9(a)  and  (b)  shows  that  the  oxygen  concentration  in  the  region 
under  the  rib  is  relatively  higher  for  the  anisotropic  GDL  than  that 
for  the  isotropic  GDL.  This  is  due  to  the  fact  that  the  large  in-plane 
effective  diffusivity  of  oxygen  in  the  anisotropic  GDL  leads  to  a 
higher  rate  of  oxygen  transport  along  the  in-plane  direction.  When 
the  GDL  is  under  inhomogeneous  deformation,  the  porosity  of  GDL 
in  the  region  under  the  rib  is  greatly  reduced,  lowering  the  effective 
diffusivity  of  oxygen  in  both  in-plane  and  through-plane  directions. 
As  a  result,  the  oxygen  concentration  is  very  low  in  the  region  under 
the  rib,  as  shown  in  Fig.  9(c).  The  variations  in  oxygen  concentration 
in  the  CL  along  the  y  direction  are  presented  in  Fig.  10.  With  tak¬ 
ing  the  GDL  deformation  into  account,  the  oxygen  concentration 
is  greatly  decreased.  Compared  with  the  case  for  an  anisotropic 
GDL  without  deformation,  the  oxygen  concentration  in  the  region 
under  the  rib  is  decreased  by  50%  when  the  effect  of  GDL  defor¬ 
mation  is  taken  into  account.  By  contrast,  a  only  10%  decrease  in 
oxygen  concentration  is  observed  in  the  region  under  the  channel, 
as  the  GDL  deformation  is  relatively  smaller  in  the  region  under 
the  channel.  Usually,  at  a  given  cell  current  density  a  lower  con¬ 
centration  of  oxygen  in  the  CL  results  in  a  higher  concentration 
polarization.  As  clearly  shown  in  Fig.  11,  the  cathode  overpotential 
is  much  higher  when  the  GDL  deformation  is  considered.  In  addi¬ 
tion,  the  profile  of  cathode  overpotential  becomes  more  uneven 
because  of  the  inhomogeneous  deformation  of  the  GDL. 

The  distributions  of  electrical  potential  at  a  cell  current  den¬ 
sity  of  1.1  A  cm-2  are  presented  in  Fig.  12.  For  an  isotropic  GDL 


Fig.  9.  Distribution  of  oxygen  concentration  in  GDL  and  CL  for  (a)  isotropic  GDL  without  deformation,  (b)  anisotropic  GDL  without  deformation,  and  (c)  anisotropic  GDL  with 
deformation. 
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Fig.  11.  Variations  in  cathode  overpotential  along  y  direction. 


without  deformation,  electron  transport  in  the  in-plane  direction 
requires  a  higher  potential  gradient  because  of  the  lower  electronic 
conductivity,  as  shown  in  Fig.  12(a).  By  contrast,  the  electronic 
potential  for  an  anisotropic  GDL  without  deformation  displays 
almost  one-dimension  behaviour  with  little  variation  in  the  in¬ 
plane  direction  due  to  the  large  in-plane  electrical  conductivity,  as 
shown  in  Fig.  12(b).  When  the  GDL  is  under  inhomogeneous  defor¬ 
mation,  much  higher  electronic  conductivities  can  be  achieved 
compared  with  the  original  value,  especially  in  the  region  under 
the  rib.  As  a  result,  electron  transport  requires  much  smaller  poten¬ 
tial  gradients  as  can  be  seen  in  Fig.  12(c).  The  electronic  potential 
drops  from  the  CL  to  the  rib  collector  are  also  plotted  in  Fig.  13. 
Clearly,  the  largest  fall  in  electronic  potential  across  the  GDL  is  only 
about  12  mV  for  the  case  with  the  GDL  deformation,  which  is  much 
smaller  than  the  values  of  60  mV  and  80  mV  for  an  anisotropic  GDL 
without  deformation  and  an  isotropic  GDL,  respectively. 

Although  GDL  deformation  results  in  an  increase  in  the  con¬ 
centration  polarization  as  the  result  of  the  increased  mass-transfer 
resistance  of  oxygen,  it  also  contributes  to  the  decrease  in  the 
ohmic  loss  due  to  the  decreased  electrical  resistance  of  the  GDL 


Fig.  12.  Distributions  of  electric  potential  in  GDL  and  CL  at  cell  current  density  of  1.1  A  cm-2  for  (a)  isotropic  GDL  without  deformation,  (b)  anisotropic  GDL  without  deformation, 
and  (c)  anisotropic  GDL  with  deformation. 
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Fig.  14.  Profiles  of  local  current  density  for  different  cathode  potentials  of  (a)  0.3  V 
and  (b)  0.55  V. 


and  the  interfacial  contact  resistance.  Cell  performance  is  finally 
determined  by  taking  into  account  both  the  adverse  effect  and 
the  favourable  effect.  Fig.  14  shows  the  current  density  profiles  at 
given  cathode  potentials  of  0.3  and  0.55  V.  Normally,  at  the  low- 
cathode  potential  of  0.3  V,  cell  performance  tends  to  be  dominated 
by  the  concentration  polarization  as  a  result  of  the  mass-transfer 


resistance  in  the  GDL.  As  shown  in  Fig.  14(a),  the  averaged  cur¬ 
rent  density  for  the  anisotropic  GDL  with  deformation  is  about  13% 
lower  than  that  for  the  anisotropic  GDL  without  deformation.  This 
is  because  severe  concentration  polarization  results  from  the  GDL 
deformation  as  discussed  earlier.  Additionally,  the  profile  of  local 
current  density  becomes  more  uneven  for  an  anisotropic  GDL  with 
deformation  due  to  the  inhomogeneous  properties  of  the  GDL.  At  a 
moderate  cathode  potential  of  0.55  V,  however,  the  ohmic  polariza¬ 
tion  gradually  becomes  dominant,  while  the  effect  of  concentration 
polarization  fades.  Although  GDL  deformation  somewhat  increases 
the  resistance  of  oxygen  transfer,  resulting  in  a  slight  increase  in  the 
concentration  polarization,  the  significant  decrease  in  the  electric 
resistance  of  the  GDL  and  interfacial  contact  resistance  greatly  low¬ 
ers  the  ohmic  loss.  Consequently,  the  overall  current  density  for  an 
anisotropic  GDL  with  deformation  shows  a  10%  increase  compared 
with  that  for  anisotropic  GDL  without  deformation,  as  shown  in 
Fig.  14(b). 

A  comparisons  of  cell  performance  for  the  three  cases  is  given 
in  Fig.  15.  Clearly,  an  anisotropic  GDL  with  deformation  results  in 
the  smallest  limiting  current  density  because  of  the  largest  mass- 
transfer  resistance  of  oxygen  in  the  GDL.  It  is  found  however,  an 
anisotropic  GDL  with  deformation  and  an  anisotropic  GDL  with¬ 
out  deformation  show  approximately  the  same  maximum  power 
density.  More  importantly,  an  anisotropic  GDL  with  deformation 
exhibits  a  better  performance  current  densities  that  range  from  0.3 


Fig.  15.  Comparison  in  cell  performance  (case  A:  isotropic  GDL  without  deforma¬ 
tion;  case  B:  anisotropic  GDL  without  deformation;  case  C:  anisotropic  GDL  with 
deformation). 


Fig.  16.  Ohmic  loss  and  concentration  loss  for  (a)  anisotropic  GDL  without  defor¬ 
mation  and  (b)  anisotropic  GDL  with  deformation. 
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to  1.1  A  cm-2.  It  is  well  accepted  that  the  cell  practically  tends  to  be 
operated  at  moderate  or  relatively  low-current  densities  in  order 
to  provide  a  relatively  high-cell  voltage  and  efficiency.  In  this  way, 
the  case  for  the  anisotropic  GDL  with  deformation  actually  shows 
better  performance  than  that  without  deformation. 

The  concentration  loss  due  to  the  mass-transfer  resistance  of 
oxygen  in  the  GDL  and  the  ohmic  loss  due  to  the  electrical  resis¬ 
tance  of  the  GDL  and  the  interfacial  contact  resistance  are  compared 
in  Fig.  16  for  an  anisotropic  GDL  without  and  with  deforma¬ 
tion.  Clearly,  at  a  given  current  density,  the  concentration  loss  is 
higher  for  the  anisotropic  GDL  with  deformation.  The  ohmic  loss 
is  virtually  negligible  from  anisotropic  GDL  with  deformation,  as 
shown  in  Fig.  16(b),  but  it  becomes  large  for  the  anisotropic  GDL 
without  deformation,  as  shown  in  Fig.  16(a).  Technically  speak¬ 
ing,  the  results  indicate  that  a  potentially  high  performance  can 
be  achieved  for  an  anisotropic  GDL  with  deformation  by  enhanc¬ 
ing  oxygen  transport  in  the  GDL.  Two  points  must  be  highlighted. 
First,  although  GDL  deformation  contributes  to  the  decrease  in  the 
resistance  of  the  GDL  and  interfacial  contact  resistance,  this  effect 
is  usually  limited.  Second,  although  GDL  deformation  leads  to  an 
increase  in  the  mass  transfer  resistance,  the  mass-transfer  of  oxy¬ 
gen  in  the  GDL  may  be  enhanced  by  using  a  new  flow-field  and 
other  techniques.  Practically,  an  optimum  deformation  should  be 
achieved  in  order  to  maintain  high  and  reliable  cell  performance. 

4.  Conclusion 

A  two-dimensional  model  for  predicting  coupled  electron  and 
two-phase  mass  transport  processes  in  an  anisotropic  GDL  of  a 
PEMFC  is  presented  in  this  work.  The  feature  of  this  model  is 
that  it  can  be  used  to  assess  the  impacts  of  the  GDL  anisotropic 
transport  properties  associated  with  carbon  fibres  and  caused  by 
GDL  deformation.  The  numerical  results  indicate  that  for  the  GDLs 
that  are  made  of  carbon  fibres  the  anisotropic  transport  proper¬ 
ties,  particularly  electrical  conductivity,  significantly  influence  the 
distributions  of  local  electrical  potential  and  local  current  den¬ 
sity,  as  well  as  the  overall  cell  performance.  The  numerical  results 
also  show  that  the  deformation  of  the  GDL  also  plays  an  impor¬ 
tant  role  in  determining  the  performance  of  a  PEMFC.  Due  to  the 
inhomogeneous  GDL  deformation,  local  transport  properties  in 
the  GDL  and  the  interfacial  contact  resistance  are  all  significantly 
changed,  which  results  in  more  uneven  distributions  of  reactant 


concentration  and  current  density.  It  is  found  that  GDL  deformation 
significantly  reduces  the  limiting  current  density  of  the  cell  because 
of  the  increased  mass-transfer  resistance  in  the  deformed  GDL, 
although  the  ohmic  loss  becomes  smaller  as  a  result  of  decreased 
contact  resistance  after  deformation.  Practically,  an  optimum  defor¬ 
mation  should  be  achieved  such  that  a  higher  limiting  current 
density  and  lower  ohmic  resistance  can  be  maintained. 
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